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Recently Hattori-Itakura have derived the full Landau-level summation form for the pho- 
ton vacuum polarization tensor in external constant magnetic fields at one-loop level. 
The tensor has three different form factors depending on the tensor direction against the 
external magnetic field. The renormalization is not trivial because these form factors are 
expressed in terms of double or triple series summation forms. We give a numerical UV 
subtraction method which can be applied to numerically evaluate the form factors in 
external constant magnetic fields. We numerically investigate the photon vacuum polar- 
ization tensor in the form of the Landau-level summation and estimate the systematic 
^H ^ errors coming from the truncation of the Landau-level summation. We find that the error 

f^^ , is practically controllable at O(10~^) level for electrons and muons in strong magnetic 

r/^ ' fields expected in heavy ion collision experiments with experimentally feasible kinematic 

parameter regions. 
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1. Introduction 

Photon vacuum polarization is a fundamental tool to access the structure of quan- 
tum vacuum. Strong external electromagnetic fields could affect the structure of 
QED vacuum and cause various nonperturbative phenomena such as pair production 
via the Schwinger mechanism, electron-positron pair production from a photon, and 
vacuum birefringence of a photon etc.^ - There have been many theoretical works to 
evaluate the vacuum polarization tensor in strong electromagnetic fields to investi- 
gate these phenomena.^ ^^ Especially heavy ion collision experiments at RHIC and 
LHC could generate very strong magnetic fields of order of eB ^ 0(m^),E2lES which 
could affect the QCD phase structure via the various chiral magnetic effects.'iS'i^ 
These effects can be experimentally looked into only if the invariant masses and the 
transverse momenta are in the accesible regions by the detectors .l^^t^il We aim at 

* ishikawa@theo.phys.sci.hiroshima-u.ac.jp 



April 15, 2013 1:29 WSPC/INSTRUCTION FILE msvl.04 HUPD-1302 



2 K.-I. Ishikawa, D. Kimura, K. Shigaki, A. Tsuji 

providing a quantitative assesment of the effects to allow evaluation of experimental 
feasibilities. 

Recently Hattori-Itakura have obtained the full Landau-level summation form 
for the photon vacuum polarization tensor in external constant magnetic fields aim- 
ing for the vacuum birefringence of a photorpSi by strong magnetic fields. Their 
form is expressed in terms of double or triple series analytically. This is based on 
the summation on the Landau-level of virtual charged particles and the form is 
very suitable to investigate the threshold structure of single photon decaying into a 
charged particle-antiparticle pair trapped in the magnetic field. The photon propa- 
gation could be a detector for the existence of very strong magnetic fields expected 
in heavy ion collision experiments at RHIC and LHCl^EllThe analytic form could 
be an important basis to discuss the polarization tensor in more realistic external 
fields>2 

In Ref . [25] they have discussed the renormalization and the UV-structure of the 
form factors contained in the vacuum polarization tensor. However their subtraction 
is not suitable to evaluate the form factors numerically because the subtraction is 
defined between the series upper limit and the UV cut-off of the subtraction integral. 
We modify their subtraction preferable to the numerical evaluation of the vacuum 
polarization tensor. 

This paper proceeds as follows. In the next section, we present the master for- 
mula for the vacuum polarization tensor in external constant magnetic fields written 
in the Schwinger's proper time integral. Our subtraction method and the Landau- 
level summation form are explained in Section [31 The convergence and the sys- 
tematic errors from the truncation of the Landau-level sumamtion are discussed in 
Section [31 The numerical results with kinematic parameters accessible in the heavy 
ion collision experiments are given in Section [SI We summarize the paper in the last 
section. 

2. Vacuum Polarization Tensor in External Constant Magnetic 
Fields 

The vacuum polarization tensor 11^"^ in external fields has the following tensor 
structure. 

Il^"'{k) = (P^"" - Pf - P^) iVo(fc) + PfNiik) + P7iV2(fc) (1) 

where k^ is the photon four-momentum and the projection tensors are defined by 

P^"^ = fc2^^- - A:^r , P|["' = fc|r?|^"-fc|'fc|(, P';^ = klT^i^^ - k^^k^. (2) 

We define the external magnetic field B is directed along z-axis and B^ = B > 0. 
The photon four momentum /c^ and the metric 77'*" are classified according to the 
direction of B as follows. 

fcff = (A:°,0,0,fc-'^) = (cj,0,0,fc,), k''^ = {0,k\k^,0) = iO,k,,ky,0), (3) 
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Tyf^^^diagll, 0,0,-1), 7?f = diag(0, -1,-1,0), (4) 

kl^{ky-{ky^cu^-kl (5) 

kl^-ikr-ikr--ikl + kl) = -kl. (6) 

We consider Dirac fermions with unit charge e > and mass m. The one-loop 
contribution to the form factors, Nj^s [j — 0, 1, 2), is given by 



TV,- = / dv dz 

47r 7_i 7o-^e 



7Vj(z,w)e"'^('^''')''-*"^(''''''^)^ - - — —e^'i 



(7) 



~ , , cos(uz) — II cot(z) sinfuz) 

No{^,v) = -^-^ —^ ^^, 

sin(z) 

Ni{z,v)^{l~v^)cot{z), (8) 

cos{vz) — cos(2) 



N2{z,v) = 2- 



sin''(z) 



, , , cosfvz) — cosfz) ,„. 

sm(z) 
1 - (1 - v^)r 
M 
where we introduce dimensionless parameters /i, r, and rj defined by 

eB kj 2(7 , fc2 

/i=^, r^—^, ?7==— , with q^-r^- (H) 

The z integration, originating from the Schwinger's proper time, should be carried 
out slightly lower line along with the real axis in the complex plane to have the 
Feynman propagator boundary condition. 

All we want to know is the value of A'j's with arbitrarily values for /z, r, q. When 
< r < 1 the both exponential factors in the integrand converge to zero on the 
lower quarter circle path with infinite radius. The integrand has no poles in the 
lower complex plane. Thus the z integral path can be continued to lower imaginary 
axis as z = —ix via the Cauchy's integral theorem. This yields the well converging 
double integration form for -/Vj's. On the other hand, when 1 < r, the z integral 
should be evaluated via the residue theorem.ESU^Zl j^ this case the line integral is the 
closed path on the quarter sector in the first quadrant of the complex plane, and 
a difficulty arises when evaluating the residue of the integrand. The poles of Nj^s 
locate on z = mr (n = 1, 2, • • • ), where the phase factor tp in Eq. ([9]) also has pole. 
This means that the residue at z = m: [n — 1, 2, • • • ) is indefinite except for the 
?7 = cases. 

When ?7 = (equivalently q — Q) case the z integration can be analytically 
performed yielding the DiGamma functions. The form for A'^i has been obtained in 
Ref. [28l Similar form can be obtained for both A^o and A''2 . 
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The Hattori-Itakura formula opens the way to evaluate iVj-'s for 77 > with 
1 < r. They used other forra factors defined by 



Xo 



-A^o, Xi = -(^i-^^o), X2^-iN2-No), 



(12) 



where the zero-field counter term is contained only in xo- For xi (X2) the subtraction 
occurs between Ni {N2) and Nq. They analytically expand the integrand of Eqs. ([5])- 
PHI in terms of C" (77) defined by 



cnv) = 



'rj 



il 



(i + n): 



{Lnv)y 



(13) 



where L^irj) is the associated Laguerre polynomial, and integrate both of the v 
and z integrals. The expansion yields the double series on n and £, where the 
indexes n and £ correspond to the Landau-level of virtual fermions trapped in the 
external magnetic field in the one-loop diagram. The zero field counter term in xq 
still remains in the integral form. The UV cancellation between the double series 
and the double integral is impossible numerically. Although xi (X2) is completely 
expanded in the double series, the location of the UV-divergence in the series is 
not aligned well between Ni {N2) and Nq. Thus the UV cancellation is not trivial 
even if they are renormalized. Therefore we need a well organized renormalization 
method suitable for the numerical evaluation. In the next section we show our 
renormalization method. 



3. Subtraction method and the Landau-level sum form 

We rearrange Eq. ([7]) to the following form. 
»i 



A^,- = / dv 



An 

a 

47r 



0~ie 
1 poo — is 

dv / dz 

-1 Jo-ie 



N. 



(z,u)e~'"''(^'^)'' ('e-"^(^^'-''^)^ - e^'i) 



iVj(z,w)e-*'^(^'")'' - - — — 1 e"'^ 



(14) 



The first line is UV finite and we can entirely expand it according to the Hattori- 
Itakura's expansion method. The second part of Nj corresponds to Nj at r = 0, for 
which z integral can be analytically continued to z = —ix yielding suitable forms 
for the numerical integration. We note that for N2 form factor the singularity of N2 
is worse than others as it contains l/sin'^(z). This yields another difficulty for the 
numerical evaluation. To tame the difficulty we replace A'^2 with 



N, = 2i- 



1 



d 



sm [z 



) dr] 



(15) 



for the first line of Eq. (fT4| before applying the Hattori-Itakura's expansion. 
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Applying the above prescription we obtain the following form for the form fac- 



tors. 



N, 



CXJ LXJ 



Att 



n 


=0 




1=0 




A^J 


ri 




/•OO 




-1 


dv 


I 


dx 



iVj(a;,v)e'^(^''')''- 



where C„ 



(2 



and 



-K^)[C;(ry)+Q"_i(77)]c;7(r,^), 

Nj{x, v) ~ —iNj{~ix, w), ipix, v) — —iip{—ix, u), 

with C\{ri) = 0. The functions, J7, G^ , Wl, and Tq, are 
F^(,r,^J)=^l[F^{r,^J)-F^{Q,^x)], 






dv 



* 



-1 
1 

1 

dv , 

-1 ^i{v;r,fj,} 

dv : r, 

-1 Sf{v;r,fj.) 



2fi 



* 



2^ 



1 



V 



S"{v] r, /i) = rv'^ — {nii)v + 1 — r + {21 + n)ji — is 



(16) 



(17) 
(18) 

(19) 

(20) 

(21) 
(22) 
(23) 

(24) 

(25) 

(26) 

(27) 
(28) 



where the z integration is performed and the Feynman's is prescription is restored 
to identify the absorptive part of these functions for the v integration. The function 
^(z) is the DiGamma function. We follow the notation for F^", G", and iJ" given 
by Ref. 1251 and the v integration can be done analytically as given in [Appendix A| 
except for 7?," . The form of 7?." is inspired from Ref. [28] in which the analytic form 
for iVi with g = has been obtained. 

The DiGamma function ^(z) has poles at z = and negative integers. When 
1 < r the argument S'"_,_i/(2m) of ^(z) in Eq. ([24]) could hit the singularities in 
integrating v. In order to extract the absorptive part of Eq. ([24]) we employ the 
recurrence formula, ^(z) = ^(z + 1) — 1/z, until the argument becomes a non-zero 
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positive number as 



1 ^ 1 

v|/(z) = *(z + 1) -- = ••• = vl,(^ + i^ + 1) _ ^ 



fe=0 



z + fc' 



(29) 



where K is a nonnegative integer chosen to satisfy z + K + 1 > 0. Thus Eq. ([M)) 
becomes 



nnr,f^) 



dv 



K 



* 



'^"+i+i<'+i(^!'''A*) 



2/i 



* 



2/x 



K 

A" 



fc>0 

-Ceihng[A,"+J 
-1 



1 
2^ 



1 - r + (2^ + n)/x - 



( |n/V(2r)| < 1 and A^^^ < 0) 
(otherwise) 



4r 



(30) 

(31) 
(32) 



The absorptive part is extracted as the sum of J-"" and the v integral can be numeri- 
cally evaluated. When \{v'^ — l)r/{2(j,)\ < 0.01 the integrand of Eq. ((30)) is evaluated 
using 8th order Taylor expansion to avoid loss of significant digits; 

*(z + dz) - *(z) ~ ¥'\z)dz + *(2)(z)i^ + . . . + *(«)(z)i^, (33) 

2 8! 

where '^^^'{z) is the poly gamma function of order j. To reduce the cost of numerical 
integrations at each £ we can use the following recurrence formula for TZ"; 



7^?(r,M)=7^^l(r,A*) + 2J■,"(r,^). 



(34) 



The form factors below threshold (r < 1) can be evaluated numerically for any 
< g by analytic continuation z = — ia; in Eq. ([7]). The form factors in these region 
have been investigated in Ref. [521 The values from the double integral are compared 
to our numerical estimates from the Landau-level summation Eq. (J16p to check the 
consistency. 

With the vanishing transverse momentum {q — 0), the z integral in Eq. ([7]) 
can be performed analytically. Karbstein et al)^ have shown the analytic form for 
Ni with g = 0, which is the integral containing the DiGamma functions similar 
to Eq. p4p . This form is valid for any r. We obtain similar analytic expressions 
for A^o and N2 with g = as given in [Appendix B| together with A^i with g = 0. 
We can check the validity of the numerical values from Eq. (|16p with g = by 
comparing to the values from the DiGamma expressions, Eqs. (|B.l|) - (jB.5|) given in 
[Appendix B[ in the case of 1 < r. Before going to the numerical evaluation, we 
discuss the convergence of the double sum of Eq. (|16|) by observing the asymptotic 
form in the next section. 
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4. Asymptotic form of the double series 

The asymptotic form for Eqs. ((25 |) - ([24l) in 1 <C ^ is given by 



When ?7 > the coefficient function C"{ri) and its derivative behave as 

^(^) e-^sin(2e,"(r?)), (39) 

ar] TTT] 

e,"(77) = 2V^-|(^n+0, (40) 

'^"=^+^: (41) 

for rj < 4k" with 1 ^ £. This is foUowed by the asymptotic form for the Laguerre 
polynomials 

L7{ri) - ^T^'-^ (K^ry)-"/^-^/''cose^(ry), (42) 

based on Bessel function expansionP^l'^ 
r^Q g and n" ^ are bounded by 

^li<Mi\-o(J^^ (forj = Oandl). (43) 

This is slowly converging series at a fixed n. For ilj £, however, it does not seem to be 
absolutely convergent since \i^2 i\ ~ '^(l/^)- The cancellation due to the oscillatory 
behavior of dCg /drj or due to the sign mixture among terms with different n could 
occur for the convergence. The worst case is that the series for ^^f, is asymptotic. 
We could not prove the convergence for ^21 with < 77 case. 

From Eq. (fT3|) the coefficient function C" [rf) and the derivative for ry = become 

C,"(0) = 5„,o, (44) 

= {l+l)5^,u (45) 






ri=0 



dC" 

—^iO)^~i2£+l)Sn,o + ie+l)SnA. (46) 

drj 
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Since these do not have damping factors for 1 ^ £ the series convergence becomes 
critical We check the convergence of the double series explicitly in the following. 
For J7q^ and 17"^ with ?/ — 0, the double sum converges as follows. 



^ C„ ^ ni, = 2 ^(J-1 - (2£ + 1)^1) 
n=0 e=o i=0 



E 

f»l 



2r 



5fiP 



^'^ 



< oo, 



£ c„ £ ni, = (j-o" -n°) + J2 2(-^£ - ^°) 



n=o e=o 



E 



_8r_ 
15m^2 



O 



P 



< oo, 



For 5I2 £ it becomes 



(47) 



(48) 



E ^" E "v = E [-4(2^ + ^)K + 8{i + ml] 



n=0 



e=o 



i=0 



«>i 



O 



< 00, 



(49) 



where the linear and logarithmic divergences are canceled among n — and n — I 
terms. Thus the double series for il^^ is not absolutely convergent and the result 
depends on the ordering of the summation. When g = (77 = 0) and r < 1 we 
observe a large discrepancy caused by the conditional convergence property between 
the Landau-level summation formula and the double integral formula numerically. 
Fortunately N2 with q ^ does not contribute to the polarization tensor as it is 
multiplied by the projection tensor Pj"' which is identical to zero. 

So far we do not discuss the convergence of the summation on n except for the 
case ?7 = 0. To check the validity of the Landau-level summation form we compare 
the value to those evaluated with the other forms numerically instead of analytically. 
The comparison is possible in the following two regions. 

(A) Double integral form in r < 1. 

(B) DiGamma form with q — 0. 

The double integral form is obtained by substituting z = —ix in Eq. ^. The 
DiGamma form is given in [Appendix B[ The integral is numerically evaluated using 
the double-exponential quadrature formula. We employ the program in Ref. [32] 
to evaluate the numerical integration not only for the above (A) and (B) but for 
Eqs. (pO]) and (fT6|) . The missing region for the validity check (A) and (B) is 1 < r 
with q ^ 0. 

Since the series coefficients C" and dCf/drj are independent from the choice of ?■ 
and only a finite set of (n, £) induces the absorptive part in J"", Q^, H" and 7?." for 
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Table 1. Parameter combinations we investigated. We use me = 0.5109989 
[MeV], m^ = 105.6583668 [MeV], m^ = 139.57018 [MeV], < fc|| < 4^ [GeV^], 
< fc5_ < 3^ [GeV^]. We take 401 (41) sample points for both r < 1 and 1 < r 
region at equal intervals. 



Case [m, eB] 



# of sample # of sample 
points for r points for q 



[a-1] 


[m^.lOm^] 


[b-1] 


[m^,ml] 


[c-1] 


[m^, (1/I0)m2 


[a-2] 


[■me,Wml] 


[b-2] 


[me, ml] 


[c-2] 


[me,(l/10)m2 



1000 




401 


31 


2000, 4000, 


8000 


41 


31 


1000 




401 


31 


2000, 4000, 


8000 


41 


31 


1000 




401 


31 


2000, 4000, 


8000 


41 


31 


1000 




401 


31 


2000, 4000, 


8000 


41 


31 


1000 




401 


31 


2000, 4000, 


8000 


41 


31 


1000 




401 


31 


2000, 4000, 


8000 


41 


31 



a finite r, the double series does not change the asymptotic form irrespective of the 
choice of r. Therefore we expect that if we have the validity in (A) r < 1 region with 
a truncated double series, the same truncated series is also valid in 1 < r region. 
The comparison of (B) provides a limited consistency check for the statement. In 
the next section we make the above comparison numerically. 

5. Numerical Results 

We employ Fortran 90 language to evaluate Eq. P^ in double precision. The double 
integral on x and v is evaluated by nesting the double-exponential quadrature for- 
mula subroutine.'^ In order to avoid the loss of significant digit in Nj{x, v) {ip{x, v)) 
near a; = we use the 9th order (11th order) Taylor expansion form for x < 0.02 
respectively. To avoid the overflow of hyperbolic functions in Nj{x,v) and ip{x,v) 
for 1 <C cc we transform them to a well organized exponential form for 10 < x. 

The coefficient functions C" and dC"^ /drj are computed using the three-term 
recurrence formula based on the Laguerre polynomials during the summation on I. 
This means that we keep several last values of L^_i and i"_2 e^c. to compute C" 
and dCl^/dri for fi"^ to avoid the full re-computation of L" at each i. The three-term 
recurrence formula and numerical method we employed is explained in [Appendix D[ 

The double series of Eq. ()16|) must be truncated at a cutoff index (rimax, ^max) for 
the numerical evaluation. The summation on £ is truncated at a ^max independent 
of n. While the summation on n is stopped when the partial sum SNj = X^fccT ^1 i 
becomes negligible compared to the current estimate of Nj provided by the following 
condition; 

{\SN^\ < 10"" and \N^\ < IQ-^'') or \AN^\/\Nj\ < 10"" (50) 

in double precision arithmetic. 
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-5 - 
-6 
-7 - 



|k, 1=0 [GeV] 

|k,|=3*Sqrt(4715)[GeV] 

ik^i=3 [GeV] 



0.0 



Fig. 1. Form factors A^o (left) and A^i (right) for muons (case [a-1]) in r < 1 with la 



1000. 



'< 0.0 



-1.0 




|k , 1=0 [GeV] 

|k,|=3*Sqrt(4715)[GeV) 

|k^|=3 [GeV] 



'< 0.0 



-1.0 





|k , 1=0 [GeV] 

|k,|=3*Sqrt(4715)[GeV) 

|k^|=3 [GeV] 



0.0 



-0.5 
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|k,|=3*Sqrt(4715)[GeV] 

|k^|=3 [GeV] 



0.0 ry 



-0.5 




|k, 1=0 [GeV] 

|k,|=3*Sqrt(4715)[GeV] 

|k^|=3 [GeV] 



200 



300 



100 



200 
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Fig. 2. Same as Fig[T]but for real (top) and imaginary (bottom) parts in 1 < r. 



We show the combination of input parameters for the form factors in Table [T] 
We choose the magnetic field strength at 0{rn^) which is expected to exist in the 
heavy ion collisions at LHC!^'^ The longitudinal and transverse momenta ranges 
we investigated are < fc| < 4^ [GeV^] and < fc5_ < 3^ [GeV^] respectively. 

Figs. [IE] show the form factors iVo, A^i and N2 with m = m^ and eB = lOmJ 
(case [a-1]). The upper limit on (. is (-max — 1000. Figs. HUH] are for electrons m — m^ 
with eB — lOrn^ and £max — 1000 (case [a-2]). Complicated threshold structures 
due to the Landau-levels are seen for 1 < r in Figs. [2]|3]and in Figs. [5]|6l The solid 
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E 
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f 
15 
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- 
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Z -1.0 



-1.5 



-2.0 






iif'T -t I 






pjVm 



|k,|=3*Sq 11(4/15) [GeV] 
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Fig. 3. Form factor N2 for muons (case [a-1]) with ^max = 1000 (top (r < 1), middle (real part 
in 1 < r) and bottom (imaginary part in 1 < r). Right panels in 1 < r arc magnification of left 
panels. 



lines in the top panels of Figs. |3] and El which correspond to A'2 with q = 0, have a 
different behavior compared to the other lines. This is because of the conditionally 
convergent property of fl^ g as explained in the last section. Thus the solid lines in 
the real part of 7V2 in 1 < r (middle and bottom left panels of Figs. -[3] and [6]) also 
contain the same systematic error. 

The truncation on the summation n is monitored as shown in Fig. [71 It requires 
360-370 terms on n for larger transverse momenta q. We also observe rimax depends 
linearly on -\/4nax resulting n^ax — 920-940 at ^niax — 8000. We observe a similar 
behavior on n„iax for electrons except for A^i. An early truncation in r < 1 for 
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Fig. 4. Form factors A^q (left) and A^i (right) in r < 1 for electrons (case [a- 2]) with £max = 1000. 
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Fig. 5. Same as Fig[4]but for real (top) and imaginary (bottom) parts in 1 < r. 



Ni is seen as it is well approximated by the lowest Landau level approximation.'231 
As we decreasing eB to m^/10, rimax increases to 5460-5500 (with £max = 1000) 
for both electrons and muons. To approach the zero field limit we must accumulate 
more contributions from higher Landau levels. Verifying the zero field limit becomes 
difficult numerically. The zero field limit for the imaginary parts with q ~ Q can be 
analytically taken as shown in [Appendix C| 

The left panels in Fig.|8]show the discrepancy between the Landau- level form and 
the double-integral form as the consistency check (A) in r < 1 for muons with eB = 
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Fig. 6. Same as Fig|4]but for N2. 



10ml and £max = 1000 (case [a-1]). The discrepancies AA'^o and ANi is at O(10~^) 
and decreases as increasing q. This is practically satisfactory level. While for AA''2, 
the discrepancy with q = (solid line) has a 0(1) error and it rapidly decreases 
to O(10~^) as increasing q. The truncation errors depend on r linearly, which is 
consistent with our asymptotic analysis. For electrons with eB — 10ml (case [a-2]) 
we observe the same behavior in r < 1 except for A^i. Since the relative truncation 
error |(AA^i)/iVi| for electrons reaches the limit of double precision accuracy, we 
cannot extract the proper r dependence for AA^i (case [a-2]). 

The right panels in Fig. |8] show the discrepancy between the Landau- level form 
and the DiGamma form in 1 < r with q = as the consistency check (B) for 
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<for;Vj(/„,,.1000) 
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Fig. 7. Jimax for inuons (case [a-1]) with £„ 



1000. 



muons with eB = lOmJ and ^max — 1000 (case [a-1]). The imaginary parts perfectly 
coincide with for all form factors. The real parts for A'o and A'^i are linearly continued 
from the left panels and still remain below 0(10^*) in the region we investigated. 
For AiV2, however, it reaches 0(1). If we extend the observation in the region r < 1 
to the region r > 1, we expect even with q > Q that AiVo and AiVi still remain 
at 0(10-4) and AiVa with \k±_\ > 3^4715 - 1.5 [GeV] remains at 0(10-^) (see 
Fig. 111). 

The q dependence of the truncation error at r = 0.8 with ^max — 1000 is shown 
in the left panels of Fig. [TU] for the case [a-1]. AA'o and AA^i behave as a linear 
function of ^^q while AiV2 behaves as a linear function of l/y^ as shown by the 
fit lines in the figures. The £max dependence of the truncation errors at r = 0.8, 
|fc_L| = 3-\/14/15 [GeV] is shown in the right panels of Fig. [10] for the case [a-1]. 
The truncation error for A^o and A^i can be fitted with (c -I- d/-\/4nax)/^max- AiV2 
can be fitted with c/V^max + d/^max- The same behavior is observed for other q in 
r < 1. This behavior cannot not be understood from the asymptotic behavior on £ 
at a fixed n because it involves the truncation effect on the n summation. 

With the global analysis for all cases shown in Table [H we find that the trunca- 
tion error can be well expressed by 



AN,/a^ [cj+d, 



ATVa/a = {c2 + d: 




(j = and 1), 



(51) 
(52) 
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Fig. 8. Comparison (A) (left panels) and (B) (right panels) for muons (case [a-1]) with £„ 
1000 



in 7' < 1 derived from check (A), and by 



ij = and 1), 



(53) 



in 1 < r derived from check (B). 

Table [2] shows the coefficients for Eqs. (f5T|) -(|53 |) obtained by fitting all data 
from the parameter sets shown in Table [TJ For electrons in strong magnetic fields 
eB — lOm^ and m^ we cannot determine ci and di properly by fitting as jAA^i/iVil 
in r < 1 reaches at double precision limit 10~^^-10~^^. We note that ANq and AA^i 
with q = give upper bound for the truncation errors since the coefficients do and 
di are negative. As seen from Table [2] cq and eg (ci and ei) are consistent except 
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Fig. 9. Figs.|8]are combined in log-log plots (case [a-1] with lyna 



1000). 



for the cases at eB = m^/10 (cases [c-1] and [c-2]). This is because the truncation 
error deviates form the function form Eq. (j53p . Thus the truncation errors for iVo 
and A'^i can be directly estimated from the comparison (B) in 1 < r region for 
sufficiently strong fields. We extend Eq. ((52|) determined by fitting in r < 1 to 
estimate the truncation error AA'^2 in 1 < r. This gives a upper bound for AA^2 
because it monotonically decreasing as increasing q. 

A practical algorithm to compute the form factors is summarized as follows; 



(1) Use double integral forms for all iVj in r < 1. 

(2) Use Landau level summation forms in 1 < r for A^o and Ni with the truncation 
control by Eq. ([55]) . 

(3) Use Landau level summation forms in 1 < r for N2 with the truncation control 
by Eq. (El 



For muons in strong magnetic fields eB = rn^-lOmJ, summation up to ^max — 
10000-20000 yields -- lO""* accuracy for TVq and Ni, and -- 10"^ accuracy for N2 
in the kinematic region; 1 [GeV] < |fc_L| and < fcn < 4^ [GeV^]. 

For form factors, especially for A'^2, in weaker magnetic fields or with more precise 
values, it becomes difficult to obtain accurate form factors with our naive summation 
method. We might need to apply series acceleration techniques. 
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Fig. 10. q dependence (left panels) and <?max (right panels) dependence of truncation errors at 
r = 0.8 (case [a-1]). 



6. Summary 

Wc investigated the vacuum polarization tensor in constant background magnetic 
fields based on the Hattori-Itakura's Landau-level summation formula with the ap- 
propriate UV subtraction method we constructed. We could reproduce the numer- 
ical values computed with the Landau-level summation form consistent with those 
with known formulae. The Landau-level summation was truncated and we estimated 
the truncation error in a range of the parameter sets for muons and electrons. In 
very strong magnetic fields eB — m^-lOm'^, we could evaluate the form factors 
with practically acceptable accuracy in the limited kinematic region 1 [GcV] < |fcj^|, 
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Tabic 2. Fit results for Eqs. |[5T1l-l53t. 



Case 


[a-1] 


[b-1] 


[c-1] 


CO 


1.820 X 10-=* 


1.794 X 10-2 


1.636 X 10-1 


do 


-2.748 X 10-4 


-7.002 X 10-3 


-1.174 X 10-1 


eo 


1.828 X 10-» 


1.877 X 10-2 


3.132 X 10-1 


ci 


2.425 X 10-3 


2.380 X 10-2 


2.116 X 10-1 


di 


-4.720 X 10"* 


-1.202 X 10-2 


-2.048 X 10-1 


ei 


2.439 X 10-3 


2.513 X 10-2 


2.652 X 10-1 


C2 


-1.324 X 10-3 


-1.231 X 10-2 


-1.203 X 10-1 


d2 


1.145 X 10-2 


3.609 X 10-2 


1.137 X 10-1 


Case 


[a-2] 


[b-2] 


[c-2] 


CO 


4.256 X 10-* 


4.196 X 10-^ 


3.830 X 10-"^ 


do 


-3.105 X 10-" 


-7.921 X 10-1" 


-1.329 X 10-** 


eo 


4.275 X 10-» 


4.391 X 10-7 


7.360 X 10-6 


ci 


- 


- 


4.959 X 10-*^ 


di 


- 


- 


-2.324 X 10-* 


ei 


5.705 X 10-» 


5.880 X 10-7 


1.076 X 10-5 


C2 


-3.096 X 10-« 


-2.879 X 10-7 


-2.820 X 10-« 


d2 


5.539 X 10-5 


1.746 X 10-* 


5.505 X 10-* 



< fcn < 4^ [GeV^] for muons and electrons. This kinematic region is accessible 
provided by a small invariant mass in the heavy ion collision experiment at RHICpH 
and LHCpSl where such a strong magnetic field exists in the early stage of the heavy 
ion collision. The propagation of a real or a virtual photon emitted in the early stage 
of the collision could receive a large asymmetry due to the direction dependent po- 
larization tensor originating from the pair creation phase space suppression due to 
the Landau-level bound states. Hadronic contributions to the vacuum polarization 
tensor must be incorporated before phenomenologically applying the propagator to 
investigate the effect of strong magnetic fields. However we expect that the polar- 
ization tensor estimated in this paper is partly applicable to prove the existence 
of strong magnetic fields via the photon propagation in the heavy ion collisions at 
LHC experiments. 
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Appendix A. Integrals for Eqs. (f25|) - (|27l) 

We follow the notations given by Ref. [25] except for the dimensionless parameters r 
and /i (correspondence to Ref. [25] is r ^^ rf,, fi ^ Br). The analytic expression for 
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Fnr,t^) = < 



In 



W 



arctan 



In 
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V 



arctan 



a-c+^/V 



2TTi 



6-2q 



(r < si") 

(si" <r < s^") 

{.sf < r) 



(A.l) 



i" = i (yTT2^±Vl + 2(^ + n)Ai)' 



(A.2) 



where a = r, b = —nji, c = 1 — r + (2^ + n)/i, and D = b^ — Aac. 
Eqs. (gni) and ^ are given by 



G7{r,^l) = -[^U^^)+n^iF^{r,^x)], 



H7[r,fi) 



En^^) ^ In 



nil 5^ — 2ac ^„ , , 

2 + -^^^^ M + -^ F\r,^i) 

2r 2a 



1 + 2^^ 



l + 2{e + n)fi 



(A.3) 
(A.4) 
(A.5) 



When evaluating these functions numerically, the naive implementation causes 
a loss of significant figures near r = 0. We use 8th order Taylor expansion forms 
when \r/{nn)\ < lO'^ forn > and |r/(l + 2^/i)| < 10"^ for n = 0. 



Appendix B. Form factors with q — 

When q = Q case we can integrate z analytically for Eq. ([7|) using the residue 
theorem and the refiection formula of the DiGamma function. The expression for 
A^i has been obtained in Ref. .28^, We give similar expressions for A'o and N2 in 
order to compare the numerical values with the Landau-level summation formula 
with (7 = as a consistency check. 
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After integrating for z in Nq we obtain 

r-l 



iVo = - 



47r 



1 



dv- 



-2w2-2(l-i;2)ln(2^) 



l + v 



1 * 



2fi 



^ J S^o{-v■,r,^,) ^\^.f Sk+^i-v■,r,f^) \ 



K 



+2E 



fc>0 



-{-2b+{c-a)bFl{r,f,) 



+ {b^ -a"- ac + a)Gl{r,^l)] + hHl{r,^i) 



(B.l) 



K 



where a = r, b = —fi, and c — 1 — r + (2fc + l)ii. The shift integer K is given by 

-Ceihng[Ai] (lAi/(2r)| < 1 and A^ < 0) 

— 1 (otherwise) 

where 5'" and A" are given by Eq. (|28l) and Eq. (|32|) respectively. 
Similarly we have 



(B.2) 



iVi 






■£U/"<'-"" 



2/i 



/< 



/C 



-M (F°(r, ^) - i7°(r, /i)) + E 2^ (F°(r, /i) - i/°(r, A*)) 

k>0 

-Ceiling[^!]] (AO < 0) 

-1 (otherwise) 



(B.3) 
(B.4) 



Nn 



a 

47r 



dv— 
-1 2 



-1- 3^2 -2(1-^2) In (2^) + 2 



^o"(^m) 



2 /^g(«^r^y ^ /^?^^7+i(^^ 

V A^ / V 2a* 

5*0 (t'; »^, m) f S^Jvyr^ -2]^ '' '^^+^ ^"' '' ^'^ 






Ai 



M 



2 * 



2/x 



2^ 



J 
'E2 



4r\ 1 
2 - y j - - 4(j + 1) + 4(j + l)Vi1Vi('^. m) 



a: 



E2 

fc>0 



4r\ 1 
2 - y j - + 2{2k + 1) - 4fc(fc + 1)//Ffci(r, /i) 



(B.5) 
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J = 

K = 



-Ceiling [A?] 
-1 

-Ceiling[Ai] 
-1 



{A\ < 0) 
(otherwise) ' 

{\^i/i2r)\ < 1 and AJ<0) 
(otherwise) 



(B.6) 

(B.7) 



Appendix C. Zero field limit of the imaginary parts of Eqs. (jB.ip 
and (|BT3l) with q = 

The vacuum polarization tensor in vacuum is written by 



a f 1 



1 



n(fc^)--|-+(^2 + -) yT7^^cot-i(yT7^^)-i 



(C.l) 



for y < I and 

a (1 



m') = ^ { 3 + (2 + - j [v/r^T7^tanh-i (yr~T7^) - 1 - ^lv^^^ 

for 1 < y with y = fc^/(2?7i)^. The imaginary part is thus 



(C.2) 



Imn(/c2) = -- (2 + l/y) Vl - l/y- 




(C.3) 



The zero field limit for the imaginary parts of Eqs. (jB.l[) and (jB.3|) can be taken as 
follows. 



lim ImA'o = lim 

eB^O eB^Q 



£'»^E 



47r 



fc>0 



lim 



{-26+(c-a)6Ffci(r,/x) 

+(62 - a^ - ac + a)Gl{r, ^i)} + bHl{r, /i) 

1-1/'- 
a TT ^ / 1 „ \ A/3 



fe>0 



/^ 



(C.4) 



where A/3 = 2/i/r and /3fc = 1 — 1/r — fcA/3. This is the rectangular approximation 
of integration and the limit leads 



lim ImA^i 

eB~¥0 



- (2 + l/r) Jl - 1/r = Imn(fcif 
6 " 



(C.5) 
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Similarly we have 



lim ImA^i = lim 



lim 

A/3->0 



47r 



-Im 



K 



/i^(2-<5o,)(FO(r,M)-if,°(r,A.)) 



fe>0 



1-1/'- 
—Kir- 



"'^EC^-Mf^lA/^ 



47r2 



fc>0 



/pk 



Air 



l-l/r 



--^.l""('-^UP 



V? 



-- (2 + 1/r) v/1 - lA - Imn(A:ff). 



(C.6) 



Appendix D. Three term recurrence for C" and dC^ /drj 

When we evaluate C" and dC™ /dr] for a large Landau level {n,m) with a naive 
implementation using the three term recurrence formula for Laguerre polynomials, 
we encounter arithmetic overflow or underflow in double precision arithmetic. In 
order to tame the numerical overflow and underflow we employ a modified recurrence 
formula with rescaling and quadruple precision arithmetic. 

We define /" and d/" satisfying the following recurrence formula; 



/o" = l, ,ri=0, dfo"-0, d/L\-0, 



(D.l) 



/; = K./ri+/S?/^2)7?, 

a^ = (2£ + n-l-77), 



/3? = {I -£-n)^ {£ -!)/{£ -1 + n), 



(D.2) 
(D.3) 
(D.4) 
(D.5) 
(D.6) 



for 1 < £. ?7 is the argument of C" and dC^/drj. /" and rf/" are proportional to 
\/£l/{{i + n)\)L"{ri) and its derivative respectively. When either of |/"| or \df^\ 
takes a value larger than 10^"° or smaller than 10^^"° during the recurrence, in- 
termediate states, (/", f('_i,f^_2, df", df^_-^,df^_-^), are rescaled by multiplying the 
inverse of max(|/"|, |(i/"|) or min(|/"|, |(i/"|) and the scaling factor is stored for 
later use below. 

The coefficients C" and dC" /di] are derived by 



dC, 



drj 



iri) 






(hndf^) - {hnfn +nign.n 



(D.7) 
(D.8) 
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„-('7+E^>llog(fc)-2Ef>T'°log(S,)-nlog(^))/2 ,^^^ ^ ^ ^. 



(D.9) 



g-(';+E^>llog(fc)-2i:f>T"=log(S,)-(n-l)log(,,))/2 ^^^^ ^ ^ _^^ 

^" " ^ g-(';+i:;!>iiog(fc)-2Ef>T'Mog(s,))/2^„_i ^^^^ 77 < 1) ' ^°'^°^ 

where S^'s are the rescaling factors stored during the recurrence. Finally the coef- 
ficients C" and dCJ^ /dt] are converted to double precision numbers. 
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